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We introduce a Lagrangian formulation of the Density Matrix Renormalization Group (DMRG) . 
We present Lagrangians which when minimised yield the optimal DMRG wavefunction in a varia- 
tional sense, both within the general matrix product ansatz, as well as within the canonical form 
of the matrix product that is constructed within the DMRG sweep algorithm. Some of the results 
obtained are similar to elementary expressions in Hartree-Fock theory, and we draw attention to 
such analogies. The Lagrangians introduced here will be useful in developing theories of analytic 
response and derivatives in the DMRG. 



I. INTRODUCTION 

The densitymatrix renormalisation group (DMRG) of 
White [l|, 0, 3 is a recent addition to the methods of 
quantum chemistry [J, [3, [y, 0, H M, [l3l ■ Unlike many 
other correlation methods the DMRG is not based on 
excitations from a Hartree-Fock reference but rather on 
a new kind of highly flexible reference function. In quan- 
tum chemistry, it has led to advances in the treatment 
of strongly interacting (i.e. multi-reference) problems. 
For molecules that are large in one spatial dimension, 
the cost of the DMRG is only quadratic in the num- 
ber of localised orbitals and it is therefore a quadratic- 
scaling multi-reference method for such systems [11| . We 
have applied a quadratic-scaling DMRG algorithm to 
study conjugated polymers [12l.ll3l|. light-harvesting pig- 
ments [14i |. and the metal-insulator transition in hydro- 
gen chains [ll| with full treatment of multi-reference cor- 
relations in as large as 100 electron, 100 orbital complete 
active spaces. 

Early formulations of the DMRG primarily used the 
language of the numerical renormalisation group that 
reflects the history of its development. Such language 
is very different from the usual language of quantum 
chemistry. However, as is now understood, the DMRG 
algorithm simply minimises the energy of a wavefunc- 
tion ansatz known as the matrix product state [3|, lid . 
llSl . Il6| . This ansatz has a very different structure from 
most quantum chemical wavefunctions and the unique 
strengths and weaknesses of the DMRG method can be 
understood from this point of view |10|] . 

In a prior publication [l0| we have presented an intro- 
duction to the DMRG from the wavefunction perspec- 
tive. In the current work, we continue along this line 
of presentation and describe simple reformulations of the 
DMRG that connect the method with well-known La- 
grangian techniques in quantum chemistry [iTI . llSl . [l9| . 
All our results are of a very elementary nature, but we 
feel there is sometimes value to writing out such things 
explicitly. In the past, Lagrangian formulations have pro- 



vided a unifled and systematic language by which to de- 
rive many results in the area of perturbation theory and 
analytic derivatives and response [ITj, il8|, ily • We believe 
that the Lagrangian formulation of the DMRG presented 
here will be useful in a similar way. 

We start in section [IT] by recalling the matrix product 
structure of the DMRG wavefunction. In section IIIII we 
write down a simple Lagrangian for the DMRG wavefunc- 
tion and the corresponding stationary equations. These 
resemble the the Fock orbital equations of Hartree-Fock 
theory and we discuss this similarity. In DMRG cal- 
culations that are based on the traditional sweep algo- 
rithm, one implicitly uses not the most general form of 
the DMRG ansatz, which contains some redundancy, but 
rather a special canonical form [lO|. The canonical form 
of the DMRG wavefunction is reviewed in the first part of 
section lTVl while in the second part we introduce the addi- 
tional constraints that have to be applied to Lagrangian 
to ensure that the minimising wavefunction is of canoni- 
cal form. We finish by demonstrating explicitly that min- 
imising the canonical DMRG Lagrangian yields exactly 
the same solution conditions on the DMRG wavefunction 
as the original sweep algorithm of the DMRG. 



II. THE DMRG WAVEFUNCTION 

Recall the full configuration interaction expansion of 
the wavefunction in Fock space. In terms of Slater deter- 
minants written in the occupation number representation 
|ni . . . Uk), where n, is the occupation of orbital i taking 
values 0, !„, 1/3, lal/3, this is 



"1": 



nin2n3---"fc 



■\nin2n3...nk), (1) 



(2) 
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In the DMRG ansatz, the expansion coefficient 
yj;nin2H3---nfc jg represented by a contracted product of 
the tensors, where each tensor is associated with the Fock 
space of a given orbital. In the context of the DMRG it 
is more usual to refer to orbitals as sites, and thus we 



refer to the tensors as site-functions. Thus we write 



with the matrix elements of the operator F[p\ defined as 



^nin2n3...nfc _ ^"i ^n2«l ^"3*2 






■^Z. (3) 



where we have used the Einstein summation convention 
that we will employ throughout this work (i.e. repeated 
upper and lower indices are contracted). The number of 
coefficients in each site-function (save for the first and 
last) is 4M^, where M is the dimension of each i index. 
(It is conventional to take the dimension of each i index 
to be the same). 

From eqn. ^ we see that the DMRG ansatz has a 
contracted matrix product structure. For this reason it 
is known as a matrix product state 0, [l3, llM 13 ■ The 
matrix product maybe used to reconstruct the Slater de- 
terminant expansion of the DMRG wavefunction (see e.g. 
|20J). The product nature is reminiscent of the orbital 
product ansatz in Hartree-Fock theory. However, there 
are some important differences. Firstly, the number of 
site-functions is the size of the basis k, rather than the 
number of electrons N. Thus the product structure of 
the DMRG is expressed in the full Fock space, not in the 
A^-particle Hilbert space. Secondly, the ansatz is a con- 
tracted product rather than a simple product. Thus cor- 
relations are introduced between the orbital Fock spaces, 
by virtue of the contraction structure of the i indices in 
the ansatz. 



III. THE DMRG LAGRANGIAN 

We can determine the best site-functions in the DMRG 
ansatz in a variational sense by minimising the energy 
subject to normalisation of the wavefunction. The corre- 
sponding Lagrangian is 



£[*] = (*|iJ|*) -£;(*|i|*) 



(4) 



To evaluate the Lagrangian explicitly with the DMRG 
ansatz we use a Fock representation of the Hamiltonian. 
The matrix elements of the Hamiltonian are written as 



(nin2 ■ ■ ■ nk\H\n[n2 . . . n').) 



TTni...nk 
n\ ...n{ 



(5) 



In terms of the site- functions, the energy term {'^\H\'9) 
becomes 

(^?'^P7''' ■ ■ ■ C ) if'V-'V (^S ^S ., ... /",-') (6) 

while the normahsation term (^|1|^) is 

Note that each term in the Lagrangian is quadratic in 
each of the site- functions. At the minimum, the deriva- 
tive of the Lagrangian with respect to the site functions 
vanishes. Then, the stationary equation satisfied by each 
site function is 



"■pV-i*p,/,"p*p-i _ fii,/,"p'p-i 



Pi'n %. 



-1 ^S 



= Ei^r^ 



F, 



'p7ipzp^i - \Vn ■■■p^ ■■■Vi^-i)^ 



H" 



= (*|iJ|*) 



(9) 
(10) 



where in the first line the struck-out symbols indicate 
that the corresponding site-functions are omitted from 
the sum, and this is denoted also by the more compact 
notation in the second line. 

The stationary equations for the site-functions are 
analogous to the orbital Fock equations of Hartree-Fock 
theory [2l| as each site-function is an eigenfunction of 
an effective site "Fock" operator F[p]. However, unlike 
in Hartree-Fock theory, the Fock operator is different for 
each site, and all site-functions possess the same eigen- 
value E. The site Fock operator F[p] may be decom- 
posed into local-site and off-site terms. Assuming the 
usual form of the electronic Hamiltonian 



H = f^a]a,+v''''^a]a\ak 



ai 



(11) 



where for simplicity we are assuming summations over 
the spin-labels of the orbitals, i.e. t^^alaj = t^'^i'^ o.la'^jcr. 
We define the local-site contribution to F[p] as 

^[pjdocai) ^ (vpliPPata, + yPPPPalalapapl^)^ (12) 
and the off-site contributions as 

^jp] (off-site) ^ (vplf Jat„^.^„^^ , = , =rt + 

v'^''^ala]akai(^not i = j = k = i = p) |*)y (13) 

This division is analogous to the division of the Fock 
operator into one-electron and two-electron Coulomb- 
exchange terms. In particular, the off-site contributions 
represent the contributions of the average "field" of all 
the sites to the local Fock operator at site p. 

IV. THE CANONICAL DMRG LAGRANGIAN 

A. Canonical form of the DMRG wavefunction 

The DMRG wavefunction as written in ([3]) possesses 
many redundant degrees of freedom. For example, given 
an arbitrary invertible matrix T, we can obtain multiple 
equivalent matrix product approximations for the wave- 
function tensor vj/"!-"'! by inserting T,T~^ in between 
two site functions, e.g. 

^ni...n, 



(8) 



"1 

n 


„;,"p»p-i,/,"p+i»p 
■ ■ ^ip ^V+i 


• ■ ■ V'"' 




"1 


^ ,np.p-l /j.Jpj.- 

V V p 


I'p'A ,;"p+iV 
ip J ^ip+i ■ 


(14) 



Thus minimisation of the DMRG Lagrangian ([4]) does not 
define the site-functions uniquely, but only up to pairs of 
transformations |10l.[l4l|. 

The original sweep algorithm used to optimise the 
DMRG wavefunction does, however, define a particu- 
lar choice of site-functions at convergence. These site- 
functions are canonical in ways which resemble the 
properties of canonical orbitals in Hartree-Fock theory. 
In Hartree-Fock theory, the canonical orbitals diago- 
nalise the Lagrange multipliers associated with orbital 
orthonormality. As we show in the next section the 
canonical site functions obtained from the DMRG sweep 
algorithm diagonalise a matrix of Lagrange multipliers 
associated with orthogonality constraints. 

Let us first recall how the sweep algorithm leads to 
a canonical form of the DMRG wavefunction and site- 
functions. We will then extract the solution conditions 
satisfied by the canonical site-functions at the conver- 
gence of the sweep algorithm. We will assume here some 
familiarity with the DMRG sweep algorithm and we re- 
fer readers to our earlier work and review for a comple- 
mentary discussion p, |lO|. (In particular, here we will 
focus on the "one-site" variant of the DMRG algorithm 

SHI,!!!). 

In the sweep algorithm the site-functions are seen as 
transformation matrices which define sets of renormalised 
bases. For example, at block configuration •i...»p_i •p 

the first p — 1 site-functions define many- 



•p+i 



• •fc 



body "left" basis functions recursively through 

Tip ^ \ ip — 2 

ni...r!p_i 

(15) 

while site-functions p + 1 . . .k define many-body "right" 
basis functions recursively through 

|rp)= Y. Kr''^^'K+^^P+^) 

+irp+i 



"p+i''p+i 



np-|_i...nfc 



(16) 



The transformation matrices are orthogonal in the sense 
that {lp-i\l'_^ = 5jf~^ and similarly for the right basis 



functions; this implies 



T ^g 

n„l 






T 



gtg-1 



'R ? 



5]f 

"-q-l 



(17) 



Using these definitions of the left and right bases {Ip-i}, 
{rp} as well as the basis of site p, {up} the total wavefunc- 
tion at the block configuration •! . . . •p-i •p 'p+i ■ ■ ■•k 



is expanded as 



l*)= E ^'-^"'''■-l^p-inprp) 

tp — \ Tip Tp 

= V L]"' . . . L]""-'^"-' X 

Z ^ ^1 tp-1 

ni...7if, 
^ip_i«prp^np+irp+i ^„,^^ |^^ ^ ^ ^ ^^^ (^g^ 

We see that the form of the wavefunction constructed in 
the sweep algorithm has a matrix product structure as 
in eqn. Q but has additional orthogonality constraints 
on the site- functions (|17p . Also, this wavefunction pro- 
vides a special meaning to the pth site-function (7'p-i"p'p ^ 
which appears as the set of expansion coefficients asso- 
ciated with the renormalised product basis {Ip-inprp}. 
We refer to the DMRG wavefunction constructed in the 
form p8)) as the site p canonical form [lO|, [lj| . 

We now recall how the the site-functions appearing 
in the site p canonical form of the wavefunction (|18|) 
are determined in the sweep algorithm. The coefficients 
(jip-^iUpTp g^j,g obtained by solving the Schrodinger equa- 
tion projected into the product basis {Ip^iUpTp} 

{l'p_,n'/p\H\lp_,nprp)C'-~'-^''-'- = EC'^-^K^'^ (19) 



These coefficients determine corresponding L, R site- 
functions at the same site p (i.e. Li' ''~^ and i?rp_i) 
as eigenvectors of appropriate reduced density matrices. 
For example, L-^^ '''^ is obtained from the eigenvectors 

of a density matrix F f ,f"^ constructed by tracing over 

^p p-i 
the right indices of the wavefunction coefficients 



-p,f tp^p_ \ j- 

n' V , Ir) 

P p—1 i' 



j^riplp-i /-ilp-\7iprpf~i 

^p'-p-i p-i P 

1 7- '^p'-p-i _ .. T nj^ljj 



wiL^ 



P Iv 



(20) 
(21) 



while Rr^ ''^^ is obtained from the eigenvectors of a den- 
sity matrix F f f obtained by tracing over the left indices 

p p 



C" 



'ip — \ Tip * ■ 



"Ci 



rp p" 



l-P-iKr' 






(22) 
(23) 



The L and R site-functions at site p do not themselves 
appear in the site p canonical form; rather we need the L 
site functions at sites 1 ... p — 1 and the R site functions 
at sites p + 1 . . . fc. But these can be obtained by solv- 
ing the effective Schrodinger equation (|19p at other block 
configurations in the sweep. Sweeping through block con- 
figurations »! 



•p-i 



»P+i ■ 



■ •fc 



for p = 1 . . . fc, and 

solving for the wavefunction coefficients C at each block 
configuration, we can obtain all the L and R site func- 
tions appearing in the site p canonical form ^8]) [a, [ifll • 
Note that any wavefunction written in the canonical 
form of one site (say p) can always be written exactly 
in the canonical form of another site (say q). In this 



sense, canonical forms at different sites are simply differ- 
ent representations of the same wavefunction [lOj. More 
precisely, given C'^'j-i"?''?^ ^^^ « '-^ at site g, we can always 
find C'p-i"''''^, i?";''""' at site p> q such that 



if' 



J- Tlqlq — 1 



.C" 



Lt) — l^lpT^p 



rjlq-iriqrq jyTipTp-x 



Jink 



(24) 



In the sweep algorithm, the conversion between the 
canonical forms of the DMRG wavefunction at neigh- 
bouring sites is known as the wavefunction transforma- 
tion 0, [13, [2J|, and it is commonly used to acceler- 
ate the convergence of the sweeps. At convergence, if 
(jip-iripTp solves the effective Schrodinger equation ([TO]) 
at site p, then the corresponding C^i-i-nqr, (determined 
through the wavefunction transformation solves the ef- 
fective Schrodinger equation ^9]) at site q. 

Let us now summarise the solution conditions satisfied 
by the site-functions appearing in the site p canonical 
form (fT5|) at the convergence of the DMRG sweep algo- 
rithm. 

1. For a specified site (p, say), the wavefunction coef- 
ficients (7'p-i"p''p+i satisfy the effective Schrodinger 
equation P^ and satisfy the normalisation condi- 



tion C^p-^'^p^p+^Ci 



p — 1 '^p'^p-\-l 



1, 



2. The L and R site- functions are each orthogonal in 
the sense of p7|) and are related to the C site- 
functions (in the corresponding canonical forms) as 
eigenvectors of the correponding density matrices 



3. The C site-functions appearing in all the canoni- 
cal forms from site 1 ... A: are related through the 
wavefunction transformation (1241). 



B. Lagrangian formulation 

Let us now show how the above conditions 1.-3. sat- 
isfied by the canonical site-functions at the convergence 
of the sweep algorithm can be obtained by minimising 
an appropriate canonical Lagrangian. We first note that 
C is constrained to have unit norm while the L, R site- 
functions are orthogonal in the sense (J17p . Thus we write 
a Lagrangian with these constraints 

£[*] = (*|7J|*) - E (C'-i«-'--Q^_,„^,^ - 1) 



q<p 






(25) 



q>p 



At the minimum, derivatives of the Lagrangian with 
respect to all L, C, R site-functions must vanish. Dif- 



ferentiating with respect to the coefficients C 



' p — 1 '^p 'j 



obtain an effective Fock eigenvalue equation for C similar 
to eqn. ^ 



"p'p-i''p 



(26) 



Comparing this with the effective Schrodinger equation 

(IT9l) we see that Fy/'^^y = {V -^n' r'\H\lp-inprp) , and 
pp— ip " f f 

thus (|26|1 is simply the same as solution condition 1 . from 
the sweep algorithm. 

Next we consider minimising L with respect to the left 
and right site-functions. In each case there are two non- 
vanishing contributions to the derivative, one from the 
energy expression l^\ll\^) and the other from the or- 
thogonality constraint. We will work out only the deriva- 
tives with respect to the left site-functions explicitly as 
similar expressions hold for derivatives with respect to 
the right site-functions. The derivative of the energy ex- 
pression is 



a/9L"''''-'(*|7j|*) 



L 



XH»^ ^C^^-inpVp^^ _^^nk \^ 

^^q I k-1 J 

-n\...n'k ( jl'^ jl'q ^ P'''*^-! 

ni...nfc \^n' ' ' ' n' i' , " ' ' ^ip-inpi-p+i • • •-"-„' 

\ i q q— I k 



H 



.Ty'^'-q^q-l^l T^q 

rinln-il' n'V 



(27) 



while the derivative of the orthogonality constraint is 



m<p 



-A\^\iq-. (28) 



and thus at the minimum, where dC/dL 









la rC 



Tlqlq^ll'g n'^lq^l ^Iq "gig-l 



(29) 



p'p, we 



Now the minimising condition ((29l) does not imme- 
diately resemble solution conditions 2. and 3. for 
the canonical site-functions from the convergence of the 
sweep algorithm. To demonstrate the equivalence, we 
first recall that any minimum of the canonical Lagrangian 
(1251) is also a minimum of the sim,ple Lagrangian ^ in 
section Hn] that did not have the additional orthogonality 
constraints. This is because we can always insert trans- 
formations as in (|14p to convert a general matrix prod- 
uct state (j3|) to a DMRG canonical form (fT8|) . and such 
transformations do not change the energy or wavefunc- 
tion normalisation appearing in the simple Lagrangian 
(|4]). Thus, given some set of i, C, R that minimise 
the canonical Lagrangian (pS)) . these all satisfy site Fock 
equations as in (|10p . Then, we can substitute the Fock 
equation (fTU|) in the energy derivative (P7)) . and we find 



or w^'ir't' 








'S'g-l'? 












(^/p-mpi-p 


,.i?-)x 






■ ■ ^v.^.K-'. ■ 


•■4i 

(30) 



Next, we transform the wavefunctions appearing in (j30p 
through the wavcfunction transformation (|24p . so that 
the C site- function is associated with site g + 1, i.e. 







(31) 



Also we observe that the C site function C^i'^i+^'^'i+^ can 
always be decomposed into the product of a matrix with 
an orthogonal matrix which we recognise as Rr^ii^''*^ 






(32) 



Finally at the minimum of the Lagrangian, all the L and 
R site-functions are orthogonal in the sense of ([TT]) and 
thus we can evaluate the contracted products of the L 
site-functions and the R site-functions appearing in (j30p 
exphcitly (substituting ^ for C"<!"5+i'-g+i) 



,( 









Di-fc-i \ _ ir'g+1 



(33) 



We can now use all these simplifications (PT]) . ([5^ . 
P3p to simplify the expression for the energy derivative 

Wy.r^l:' dSOl). We find 



the density matrix in the sweep algorithm TJ^f "" (con- 



structed from C* '"I"''''') in a simple way using eqn. ([32 






(37) 



Next, we multiply ((361) ^ii both sides with the L site 
function 






/-"g'g-i 



L 






(38) 



Substituting in ((57|) for the first bracketed term and us- 
ing the orthogonality of the L site-functions p7)) for the 
second bracketed term this becomes 



77ip"g'g-l /-"g'lj-l _ j-nglg- 



>£' 



(39) 



This is now identical to the density matrix eigenvector 
condition (P5)) (up to a multiplicative factor of E) if we 
simply perform a unitary transformation to diagonalise 
fi such that 



Mi 



^w'^fJ 



;».?'? 



(40) 



Thus we have arrived at our final result: the density ma- 
trix eigenvector condition of the sweep algorithm, which 
defines the L and R site functions in the canonical form 
of the DMRG wavcfunction, is equivalent to minimising 
the canonical Lagrangian (pS)) . up to a unitary transfor- 
mation of each site function which does not affect the 
energy but which diagonalises the matrix of Lagrange 
multipliers /i. The DMRG site-functions are thus "canon- 
ical" site functions in a way analogous to the canonical 
Hartree-Fock orbitals, which are similarly obtained from 
any energy minimising set of orbitals, by performing a 
unitary transformation to diagonalise the orthonormal- 
ity constraints |2ll |. 



^^n'V A' — "n'l' ,^l' 

gq-lg qq-1 <i 



(34) 



-.lo 



where D^J is a density matrix built from the coefficients 






(35) 



and the Lagrangian minimising condition from the energy 
derivative (^5]) becomes 



r)'g r g _ J-q t g 



(36) 



This minimising condition (|36p is in fact a density matrix 
solution condition very similar to conditions 2., 3. aris- 
ing from convergence of the DMRG sweep algorithm. To 
see the connection explicitly we recognise that the den- 
sity matrix Dj? (constructed from c'«''«+i) is related to 



V. CONCLUSIONS 

We have shown that we can write down Lagrangians 
which on minimisation yield optimal density matrix 
renormalisation group (DMRG) wavefunctions in a vari- 
ational sense. In particular, we have demonstrated the 
equivalence between minimising a canonical form of La- 
grangian, and converging the DMRG energy through the 
original sweep algorithm, up to certain unitary trans- 
formations of the variational parameters in the DMRG 
wavcfunction which leave the Lagrangian invariant. With 
an increasing understanding of the DMRG from a wave- 
function ansatz perspective, it is natural to look towards 
developing analytic derivative and response techniques 
as are available for other kinds of wavcfunction ansatz 
in quantum chemistry. The results presented here are a 
first step in that direction. 
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